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A system of two coupled ensembles of phase oscillators can follow different routes to inter-ensemble 
synchronization. Following a short report of our preliminary results [Phys. Rev. E. 78, 025201(R) 
(2008)], we present a more detailed study of the effects of coupling, noise and phase asymmetries in 
coupled phase oscillator ensembles. We identify five distinct synchronization regions, and new routes 
to synchronization that are characteristic of the coupling asymmetry. We show that noise asymmetry 
induces effects similar to that of coupling asymmetry when the latter is absent. We also find that 
phase asymmetry controls the probability of occurrence of particular routes to synchronization. Our 
results suggest that asymmetry plays a crucial role in controlling synchronization within and between 
oscillator ensembles, and hence that its consideration is vital for modeling real life problems. 

PACS numbers: 05.45.Xt, 89.75.Fb, 87. 19. La 



o 

< 



> 
m 

(N 

^' 

O 

o 



I. INTRODUCTION 

Ensembles of coupled oscillators are ubiquitous in na- 
ture. They arise in diverse areas of science including 
physics, biology, chemistry, neuroscience, social, elec- 
trical and ecological systems. Examples include syn- 
chronous emission of light pulses by populations of fire- 
flies [l|, synchronized firing of cardiac pacemaker cells 
01 , synchronization in ensembles of electrochemical oscil- 
lators [1, IJi] , both short- and long-range synchronization 
in the brain (within and between neuronal ensembles) 
1) B 13' emission of chirps by a population of crickets 
8|, and synchronous clapping of audiences in auditoria. 
Research into the dynamical properties of large ensem- 
bles of this kind has been a subject of intense interest 
since the 1960s H, [H [H, [ll. Mean-field theory 
facilitates the study of such ensembles by reducing the 
dynamics of a number of oscillators to the dynamics of 
their mean field, i.e. effectively of a single oscillator. In 
principle, each oscillator in the ensemble contributes to 
the dynamics of the mean field, so that the collective dy- 
namics of the entire ensemble can be represented by the 
dynamics of the mean field. This approach has a good 
analytical background that enables identification of bifur- 
cation boundaries and stability criteria for understanding 
the synchronization dynamics of the ensemble. Although 
the mean field approach suggests consideration of the dy- 
namics of just one oscillator in place of the ensemble dy- 
namics, recent research has identified new phenomena 
such as intra-ensemble and inter-ensemble clustering 14 1 
that can only be understood in terms of ensembles. Thus 
one should expect to model natural systems comprised 
of interacting entities as ensembles of coupled oscillators, 
rather than always approximating them as a single oscil- 
lator. 

Synchronization, or concurrence between oscillatory 
systems, is a remarkable phenomenon that is often in- 



escapable for coupled oscillators. Phase synchroniza- 
tion was first reported by the Dutch physicist Chisti- 
aan Huygens well back in the 17th century based on 
his observation of two pendulum clocks that persisted 
in precise antiphase, seemingly indefinitely. Thereafter, 
the phenomenon of synchronization has been studied 
theoretically JllJTlfll [Tel, [l7j and experimentally 
S i, i, i, 0, [M Ea, EOI in great detail. It is well 
known that the control of synchronization in natural sys- 
tems dll, [2^ [iSl is of great important. The occurrence 
of synchronization is very irnportant for e.g. lasers and 
Josephson- Junction arrays [3, |2^, cardio-respiratory 
synchronization [2^, Wk or temporal coding and cogni- 
tion via brain waves [H, [H, [sO, HH, HI]. However the 
emergence of synchronized oscillations can also give rise 
to undesirable effects, as in the case of epileptic seizures 
(33I [34 II , Parkinson's tremor [H, [3^ , or pedestrians on 
the Millennium Bridge p^ . 

In real systems, the interactions between the oscil- 
lators are often asymmetric. Examples include cardio- 
respiratory ^37., 38] and cardio-^ (EEG) interactions 391, 
interactions among activator- inhibitor systems y, [4^, 5l|, 
[42} . coupled circadian oscillators and the interac- 
tions between ensembles of oscillators in neuronal dynam- 
ics [H,!!!, li^. Neglecting coupling asymmetry, i.e. as- 
suming symmetric interactions, is an approximation that 
may simplify the analysis but which may also lead to a 
model that fails to describe important phenomena occur- 
ring in the system. We have already reported [14] novel 
global clustering phenomena, and novel routes to inter- 
ensemble synchronization that occur only in the case of 
asymmetrically interacting systems. It is evident, there- 
fore, that explicit consideration of asymmetry in the in- 
teraction may be essential to create a realistic model. 

In this paper, we supplement the preliminary account 
[Hj of our investigations of two asymmetrically interact- 
ing ensembles of oscillators by providing additional detail 
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of the different synchronization regimes, and we extend 
it by reporting the effects induced by noise asymmetry. 
We thereby emphasize the importance of asymmetry - in 
couphng, noise and phase - in such systems. We show 
that it is the couphng and phase asymmetries that con- 
trol their synchronization. We also report the occurrence 
of certain novel routes to inter-ensemble synchronization. 
We show that these routes are characteristic of asymmet- 
rically interacting ensembles of oscillators and that they 
cannot occur in systems where the interactions are sym- 
metrical. These results yield new insights into how syn- 
chronization arises in coupled oscillator ensembles. This 
understanding is an essential prerequisite for the devel- 
opment of control schemes, paving the way to possible 
ways of controlling synchronization in real systems. 

We introduce the model of asymmetrically interacting 
ensembles of oscillators, and define their mean field, in 
Sec. ini In Sec. IIIII we discuss analytically the stability 
of the incoherent (i.e. unsynchronized) state in the ther- 
modynamic limit and consider how it can be modelled 
numerically. Sec. IIVI defines the five distinct synchro- 
nization regimes that we have identified, and discusses 
in turn how each of them is influenced by asymmetry in 
coupling, noise, and phase. The several routes followed 
to synchronization, and between different synchroniza- 
tion regimes, are discussed in Sec. [V] Finally, in Sec. |Vl] 
we summarize the main results and draw conclusions. 



II. COUPLED PHASE OSCILLATOR 
ENSEMBLES 

The energy emitted or absorbed by an individual os- 
cillator in the ensemble will alter the physical states of 
the neighbors to which it is coupled; in particular, the 
periods of its neighbors are altered (either lengthened or 
shortened). The way in which the period is altered de- 
pends on the state of the neighbouring oscillator at the 
moment when it receives the impulse. One of the com- 
monest scenarios to consider is an ensemble of nonlinear 
oscillators evolving in a globally attracting limit cycle of 
constant amplitude. Such oscillators are called limit cy- 
cle or phase oscillators. If they are coupled in such a way 
that they will not be perturbed sufficiently to leave their 
limit cycles, then one degree of freedom is enough to de- 
scribe the system dynamics. Let us consider a system 
of two asymmetrically interacting ensembles of oscilla- 
tors (AIEOs). Their phase dynamical equations can be 
written as p^] 



.(i,2)_^^n,2)_^ 



Ar(i.2) 



B 



Ar(2,i) 



3) E M#''^-^r + «(3))-^,p)(0. (1) 



within (intra-), and between (inter-), the ensembles; / 
and g are 27r-periodic functions that describe coupling in 
the ensembles. The fact that ^ A^"^^ implies that 
the oscillators in the ensembles are asymmetrically cou- 
pled. 9\ ' are the phases of the zth oscillator in each 
ensemble and TV^^'^) refer to the ensemble sizes; we take 
= Ar(2) = jsf. From Eq. ([T]), it is obvious that each 
oscillator will run at its own characteristic frequency uji 
when uncoupled. However when coupled, there tends to 
arise a collective behavior in the ensemble. Depending 
upon the strength of the coupling parameters, the os- 
cillators either partially or completely synchronize. The 
emergence of synchronization is spontaneous beyond a 
critical value of the coupling parameter. 

(1 2) 

The rj] ' are independent Gaussian white noises with 
(Tyf '')(<)) = and {vl'-'^\t) r;f''^'(t)) - 2i^(i'2)^(i „ 

t')Sij and X^^'^^ are the noise intensities; K'^^^ ^ K'-^^ 
represents noise asymmetry. Phase asymmetry is intro- 
duced by phase shifts < q;(1'2,3) ^ ^/2. The primary ef- 
fect of the phase asymmetry is to synchronize the oscilla- 
tors to an entrainment frequency that differs from a sim- 
ple average of their natural frequencies. Such asymmetry 
is widespread in natural systems like heart cells and 
the cardiorespiratory interactions [13, HI] . Phase asym- 
metry is used to model synaptic information and time 
delays in neuronal networks and also in the phase reduc- 
tion of nonisochronous oscillators [3]. The natural os- 
cillator frequencies uj^^''^"^ are assumed to be Lorentzianly 
distributed as g^^'^^uj) = ^(7^ + (^(i-^) - Cj(^-^^yf)-^ 
with central frequencies w^'^, and 7 is the half-width at 
half- maximum. 



A. The Mean Field 

When iV ^ 00 in the thermodynamic limit, each oscil- 
lator in the ensemble can be regarded as being coupled 
to the mean field. Thus for infinitely many oscillators, 
synchronization can conveniently be defined and charac- 
terized by a mean-field (order) parameter 



^(l,2)g*^<^-=' 



1 ^ 

N ^ 



The interactions are characterized by coupling parame- 
ters A^-^'^^ and B to quantify respectively the interactions 



Here 'ip^^'^^{t) are the average phases of the oscillators 
in the respective ensembles and A^'^^{t) provide mea- 
sures of the coherence of each oscillator ensemble, which 
varies from to 1. The amplitude of each order pa- 
rameters r'^'^) vanishes when the oscillators in the corre- 
sponding ensemble fall out of synchronization with each 
other, and is positive for synchronized states, thus char- 
acterizing intra-ensemble synchronization. When Sip — 
^(1) —ipi^) ~ constant the ensembles are mutually locked 
in phase, defining the state of inter-ensemble synchro- 
nization. Geometrically, if we consider the phases of all 
the oscillators to be moving on the unit circle, then the 
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mean field is the centroid of all the phases. With this 
characterization, we show that an increase of the cou- 
pling strength between two ensembles that are synchro- 
nized separately does not immediately result in their mu- 
tual phase-locking. Rather, phase-locking occurs through 
either one of two different routes: in Routc-I the oscilla- 
tors in the two ensembles combine and form clusters; in 
Route-II one of the ensembles desynchronizes while the 
other remains synchronized. Further, there also exists 
the possibility that phase-locking between the ensembles 
cannot occur at all. 



III. STABILITY OF THE INCOHERENT STATE 
IN THE THERMODYNAMIC LIMIT 

In the limit — > cxd, a density function can be defined 
as p^^'^^ {9,t,uj)diL!d9, to describe the number of oscilla- 
tors with natural frequencies within [u, uj ~\- duj] and with 
phases within [6*, 9 + d9] at time t. For fixed uj the distri- 
bution p*^^'^^ (6*, i, cj) obeys the evolution equation 

dt 09^^ ' 89^ ' 



where v^^'^^ are given by 



,(1.2) ^ ^(1,2) _ ^(1,2) / / ^(1,2)(^) 

Jo J ~oo 



2Tr 



X f{9~(j) + a'^^'^^)p^^^'^\(j),t,uj)duj-B I d9 

Jo 

poo 

X / g^^'^\u)h{9 - (t) + a'^^'^)p^^'^\4),t,oj)du. 



The function p^^'^^ {9,1,0;) is real and 27r periodic in 9, so 
it can be expressed as a Fourier series in 9 

oo 



= ^+p['''^e''' +c.c + 7ji9,t,u;), 

2,71 

where c.c is the complex conjugate of the preceding term 
and r]{9,t,uj) denotes the 2nd and higher harmonics. 
Substituting p'-i'^) [9, t, lo) into the evolution equation, we 
get 

+ (^^i;(l'2)+^2^(l,2))^(l,2) 

oo 

= 2*;^^(a..p|lf (2) 

where pLf ^ = P*i^^'''\ tZ;(i'2) ^ ^(1,2) _ (^(1,2)/^ + Bh^) 
and a, = (A(i.2)^.fc«u.^)^^^^(i,2)^ _^ 5^»<3)^^^^(24)^)^ 
The linearized form of Eq. ^ reads as 



pk 



= -(jA:cI;(i'2)+fc2^(i,2))^(i,2)^^^^^^ 



where the Fourier components for |Z| > fc are neglected 
since I — ±fc are the only nontrivial unstable modes, 
po = 1/2-K is the trivial solution corresponding to inco- 
herence, and /fc and are coefficient of the Fourier series 
of functions / and h. Here (•) represents the average over 
the frequencies aj^^'^) weighted by the Lorentzian distri- 
bution g(i'2)(w). Solving Eq. ^ we get 



(1,2) 7,(1.2)/ X 

Pk =K (^) 



o{\p\). 



(4) 



Substituting the above equation back into Eq. ([3]) we 
find 



^(1.2) 



(C) 



{A(^^'Hb^k'"^)+B{b^^''^)) 

(Afe -)-ifcw(l'2) ^_fc2^(l,2))' 



(5) 



where A(i'2) = ikA^^''^) fke'''°'^^'^\ B = ikEhkc' 
The integrals in this equation can be written as constants 
are to be determined in a self-consistent 
manner. Thus, for the assumption 



C^l,2) 



b^^^'\u')g'^''^\j)du' 



Eq. ([5|) for fe^^'^^ becomes 



1.(1.2)/ N 

^k = 



(^(1.2)^(1,2)^^^(2,1)) 

(Afe -f ifcw(i'2) +/^2J^(l,2))■ 



(6) 



(7) 



This on substitution back into Eq. (O results in the 
following characteristic equation 

1 = ^(i^mi + A(2)m2 - (^(1^(2) -52)mim2, (8) 

where m, = j^^{g(^\oj)duj) / {Xj, + -I- Pk'^'^), i = 

1,2. The eigenvalues obtained from ([8]) are 

1(1) +1(2) 1/ 

Afe± - ^ j-k^K-ikLu±U4B^ + AA^ 

- k{Auj - ikAK){ALu - ikAK + 2iAA)Y , (9) 

where K = ^^'"+^'") , AK = ifd) -if(2), Ac. = c^d) - 
cD(2), AA = (A(i) - A(2)) Q = (w(i) + cD(2))/2. 

For a detailed analysis of the above equation, we 
specify sinusoidal forms for the functions / and h as 
{f,h}{9) ~ s'md. Therefore the eigenvalue equation ^ 
becomes 

A± = -j-K+'^e'"'±^{^e^'°' -A{AK + iAuj)e'°' 



-Auj^ + AK^ + 2iAujAK)^ 
or equivalently we have 



A± 



K 



7- 



fe*" ± i(p2 +g2)igiiC 



lUJ, 



(10) 



p>0 



-K --f- 



g2-jlgjlC 



(11) 



lUl, 



p<0 
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FIG. 1: Theoretical B-^ bifurcation diagram for a = Q, Atj = 
= 0, with (a) = = 1, and (b) yl'i' = 1.2, 
j4'^' = 1.0. The different synchronization regimes are as fol- 
lows: NS, no synchronization; S, synchronization with single 
entrainment frequency (reached via a single Hopf bifurcation) ; 
D, synchronization with two entrainment frequencies (two 
Hopf bifurcations); SI, both the ensembles entrained to a sin- 
gle frequency; S2, either of one of the ensembles synchronized 
with single entrainment frequency; Dl, the two ensembles be- 
have as one, with the oscillators in each ensemble entrained to 
either of the two distinct frequencies; D2, synchronization in 
both the ensembles separately with two entrainment frequen- 
cies. Regardless of symmetry, the notations S and D represent 
respectively single or double frequencies occurring after one 
or two Hopf bifurcations. The boundary between regimes NS 
and S/D represents 7c+. The lines of *s represent the nu- 
merically determined bifurcation boundaries for r > 0.7. For 
comparison, a numerical boundary (•s) for r > 0.5 is also 
plotted in (a). Insets show the frequency distributions (also 
obtained numerically) for the indicated regions; their ordi- 
nate axes represent oscillator counts in thousands. Note that 
the occurrence of perfect synchronization with (a) 2000 and 
(b) 1000 oscillator groups cannot be expected throughout the 
whole of each indicated region. The line I-II in (b) is one of 
the routes to synchronization discussed in Sec. FVl 



where = a, z = 1,2.3, k = ^(^^ + A^^), A = 
(Ad) - A(2)), ^ = + ^2)^ ^ ^ tan-i(§), P = 

^cos(2a) + 74[Aa;sina - A/^cosa] - Atj^ + AK"^ , q = 
^ sin(2Q;) - AAlu cos a - AAK sin a + 2AujAK. The re- 
sultant bifurcation diagram is shown in Fig. [TJ It is dis- 
cussed in detail below, in Sees. IIV Al and IIV CI Here it is 
obvious that for the case when phase asymmetry is ab- 
sent, when A^^) = B = 0, the characteristic equation ([8]) 
reduces to the characteristic equation of the Kuramoto 
model derived by Strogatz et. al. [lol.[lll|. 

A. Numerical considerations 

To investigate the system numerically, we use a Runge- 
Kutta fourth order (RK4) routine for solving the model 
equations with a time step of 0.01 (we have confirmed 
that the results are not affected by decreasing the time 
step below 0.01). We take N = 1000 in each ensem- 
ble and the initial phases of the oscillators are assumed 
to be equally distributed within the interval [0,27r]. As 
a signature of synchronization, we take the condition 
Re{X±) > in the case of the analytic treatment. For 
the numerical experiment, we set r*^^'^) > 0.7 for intra- 
ensemble synchronization in the corresponding ensem- 
bles, and a constant Sij^ for inter-ensemble synchroniza- 
tion as the conditions. The numerical condition for intra- 
ensemble synchronization, that r*-^'^) > 0.7, may at first 
seem too strict when compared with the analytic condi- 
tion that A^'^') > 0. However, there are certain differ- 
ences between analytic and numeric considerations that 
make this choice reasonable. Mainly, N is finite for the 
numerical experiment, whereas analytic conditions are 
derived in the limit N — > oo. Further, the analytic 
and numeric bifurcation boundaries (discussed later) are 
found to match quite closely for this choice of the nu- 
meric threshold for r^^'^\ We have plotted the numerical 
boundary for r — 0.5 along with r = 0.7 in Fig. [TJa) to 
illustrate this. 



IV. SYNCHRONIZATION REGIMES 

We have identified analytically the possibility of five 
distinct dynamical regimes [l^ : 

• NS: the region of no synchronization or incoherence 
(steady state). 

• SI: the region of global (inter-ensemble) synchro- 
nization, in which the oscillators of both ensembles 
are all entrained to the same frequency. 

• S2: the region where there is synchronization 
within one ensemble but not the other. 

• D2: the region of synchronization within both en- 
sembles, separately and independently, with two 
different entrainment frequencies. 
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• Dl: a global regime in which the two ensembles 
behave as one, but the oscillators within each en- 
semble are entrained at either one of two distinct 
entrainment frequencies. We will call this phe- 
nomenon inter- ensemble clustering. 

Regions S2 and Dl cannot occur when coupling and noise 
asymmetries are absent [H, (see Fig. [T](a)). In the 
following subsections, we will consider how these syn- 
chronization regimes are affected by coupling, noise and 
phase asymmetries, respectively. 



A. The effect of coupling asymmetry 

Consider Fig. [Ijb) for the case a ~ and K'-^'> ~ 
K^'^'^ = 0, when ^ — Aw^ > 0. If we start from the state of 
no synchronization (region NS), and decrease 7 for fixed 
B > 1, the incoherent (steady) state becomes unstable 
via a single Hopf bifurcation. Thus the system enters into 
the region SI from NS (crossing 7c+) and the ensembles 
entrain to a single frequency f2+. With further decrease 
of 7 below the 7c_ line in the region Dl (in Fig. [1]), a new 
entrainment frequency emerges through a second Hopf 
bifurcation. In this region, the oscillators from the two 
ensembles combine and form two clusters (inter-ensemble 
clustering) oscillating with two frequencies. 



n± = -Im(A±) = ±(1/2) {i ^ AiJ f + Auj 



1 

X sin I - tan 
2 



AAw/(^- Aw^) +tj.(12 



The lines 7c± in Fig. [T] are obtained by imposing the 
condition i?e(A±) in Eq. (fTT|) . 

Thus in this region the order parameters r^^'^^ either 
fluctuate in a quasi-periodic manner or have complicated 
dynamics (see Figs. [SJa) and (c)). This is because each 
ensemble has two clusters oscillating with different fre- 
quencies (see Figs. Ufa) andO. Thus, the behavior of 
the order parameters r^^'^^ is quite subtle. Since r^^-^'^ 
measure only the amount of synchronization within an 
ensemble, a decrease in r^^'^-* will not necessarily mean 
desynchronization. Rather, for a sufficient value of the 
coupling parameters, a decrease in r^^'^^ represents a sig- 
nature of inter-ensemble synchronization: clustering cor- 
responds to the occurrence of desynchronization within 
an ensemble because some of its oscillators tend to syn- 
chronize with the other ensemble. 

Again looking at Fig. [Hb), when B < 1, a decrease 
in 7 takes the system from region NS to region S2 by 
crossing the line 7c+ through a single Hopf bifurcation. 
Further decrease in 7 causes the system to cross the line 
7c- and, via another Hopf bifurcation, enter into region 
D2 where there are two entrainment frequencies f2± . The 
latter can be calculated from Eq. I|lip . In region S2, 
intra-ensemble synchronization can occur in either one of 
the ensembles, depending upon whether A^^^ or A'^'^^ is 
greater; in Fig. [^b), since A'^-* > A'^'^\ synchronization 
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FIG. 2: Time variations of the coherence parameters r^^^ 
(grey), r'^' (black) and the phase difference (St/j, as obtained 
from numerical simulations. Parameter values are: (a), (b) 
B = 1, a = 0.23 and (c), (d) B = 1, a = 0.47 corresponding 
to regions Dl and D2 respectively of Fig. ^ih) as traveling 
along the line I — II. Note that in region Dl the order pa- 
rameters display no synchronization. 
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FIG. 3: The distribution of frequencies in region Dl for the 
same parameter values as in Fig. [2Ia): (a) first ensemble; (b) 
second ensemble. The splitting of the first frequency compo- 
nent into two almost indistinguishable sub-components cor- 
responds to a discrepancy between numerics and analytics, 
attributable to approximations (see text) made in the former. 



occurs in the first ensemble with the second ensemble 
remaining incoherent. Note that, on increasing B (for 
fixed 7) while in region S2, the condition ^ — Atj^ < 
is violated and the ensembles enter into the phase-locked 
region SI. In region D2, the ensembles synchronize sepa- 
rately to two the locking frequencies (unlike region Dl 
where the ensembles combine and synchronize to two 
locking frequencies given by Eq. ((12])). 

The corresponding [B — 7) bifurcation diagram for the 
case A^^'>=A^'^'^ is plotted in Fig. [T](a) to show the differ- 
ence between these two cases. Region D represents intra- 
ensemble synchronization which occurs through a degen- 
erate Hopf bifurcation (similar to the route to D2) with 
entrainment frequencies Vl± = =f(1/2)(Aw^ — B^)^ + Cj 
and S represents inter-ensemble synchronization through 
a single Hopf bifurcation (similar to the route to SI) with 
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same frequency fl = lu. Note that regions S2 and Dl 
cannot arise for the symmetric coupUng case and that 
these two synchronization regimes are therefore induced 
by couphng asymmetry. 

The presence of two entrainment frequencies in region 
Dl can be seen by fooking at the frequencies into which 
all the individual oscillators are grouped as shown in Fig. 
[3] (since the order parameters do not reveal this synchro- 
nization phenomenon). The inter-ensemble clustering 
that occurs in this case is quite different from the forma- 
tion of clusters in a single ensemble 10, 48, 49] - here the 
oscillators in two different ensembles combine and form 
clusters. The occurrence of this phenomenon provides a 
new insight into possible ways of controlling synchroniza- 
tion in more realistic situations (considering asymmetry) 
like neural networks where some neurons from one ensem- 
ble (say cortex) tend to synchronize with other ensemble 
(say thalamus) creating desirable (temporal coding) or 
undesirable effects (as in the case of epileptic seizures). 
For instance, in a thalamocortical model of the neuronal 



synchronization mechanisms during anaesthesia [32l |. we 
found that the transition from deep to light anesthetized 
state occurs as a result of a fraction of the thalamic neu- 
rons entering into synchronization with the cortex, at the 
same time losing synchronization within its own ensem- 
ble. The clustering that occurs in this case is desirable 
in the sense that it favours coding of sensory information 
and helps the brain to resist the effects of anaesthesia 
and successfully maintain consciousness and cognition. 
Without coupling asymmetry, these phenomena would 
not occur. 
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FIG. 4: Theoretical B-^ bifurcation diagram for asymmetric 
noise with A*^' = A'-^^ ^ 1, AK = 0.25 with K'^^ = 0.3 and 
A''^' — 0.05, Au) — 1. Note that the synchronization regimes 
S2 and Dl emerge in the presence of noise asymmetry even 
for symmetric coupling (cf. Fig. lUa)). 




B. The effect of noise asymmetry 

It is well known that real physical systems are in gen- 
eral subject to noise. Here, we regard as "noise" any 
kind of random fluctuation in the system, whether orig- 
inating internally or externally. Synchronization effects, 
induced and modified noise, are one of particular interest 
[53,IIll[I3,[5l,[5l. 

When asymmetric noise in introduced into a system 
with asymmetric coupling, the bifurcation regimes re- 
main the same in the presence of coupling and phase 
asymmetries. There may be changes in the boundaries of 
the respective regions and their entrainment frequencies. 
However, for the case of symmetric coupling, asymmet- 
ric noise can induce the phenomenon of global clustering. 
We have already seen in Sec. IIV Al that inter-ensemble 
clustering phenomena (region Dl) cannot occur in a sys- 
tem with symmetric coupling and symmetric noise. Fig. 
[5] plots the individual oscillator phases, determined nu- 
merically, indicating the transition from SI to Dl induced 
by noise asymmetry. When A^^^^ — A'-^^ — 1.4, B — 1, 
7 = 0.05 and AK — the system is in region S (cor- 
responding to Fig. [1] (a)) where synchronization occurs 
in both the ensembles with one entrainment frequency. 
This can be seen from the top panel of Fig. O where os- 




FIG. 5: Time evolution of the oscillator phases in the first 
(grey) and second (black) ensembles. Parameter values are 
Al ^ A2 = 1.4, B = 1, 7 = 0.05, and either AK = 0, 
with iff^) = 0.2 and A'*^' = 0.2 (top) or AK = 0.25, with 
K^^^ = 0.2 and A''^^ = -0.05 (bottom). Thus the top and 
the bottom panels represent respectively the synchronization 
regions S and Dl, induced by noise asymmetry. 



cillators from both ensembles lock to form a single major 
cluster. On the other hand, when AK — 0.25, the com- 
bined system of the two ensembles synchronize to two 
main clusters, each of which comprises a fraction of the 
oscillators from both ensembles (see Fig. [S] (bottom)), 
representing region Dl. Thus it is becomes obvious that 
asymmetric noise can in some ways imitate the effects of 
asymmetric coupling when the latter is absent. 

The S2 region also appears in this case, induced by 
noise asymmetry. Here too, depending upon whether 



7 



(a) 



(b) 



2 

2 <»»> i^iiiiyv 



400 



450 



500 400 



450 



500 



(c) 



(d) 



400 



450 



500 400 
time 



450 



500 



FIG. 6: Noise asymmetry-induced synchronization regimes 
obtained numerically for the case of symmetric coupling. 
Black and grey lines represent the time evolution of the os- 
cillator phases in the first and second ensembles respectively, 
for the same parameter values as in Fig. |3]and: (a) B — 1, 
7 = 0.1; (b) B = 0.5, 7 = 0.1; (c) B = 1.2, 7 = 0.4; and (d) 
B = 0.6, 7 — 0.4. Panels (a)-(d) represent the synchroniza- 
tion regimes Dl, D2, SI, S2 respectively. 



AK is positive or negative, synchronization occurs either 
in the second or the first ensemble respectively, similar to 
the case when region S2 arises in the presence of coupling 
asymmetry. Thus noise asymmetry plays a similar role 
to coupling asymmetry for the symmetric coupling case, 
and the {B — 7) bifurcation diagrams [ijb) and [4] look 
similar. Fig. [S] depicts the results of numerical investi- 
gation of all the synchronization regimes in the presence 
of noise asymmetry corresponding to the analytical bi- 
furcation diagram in Fig. |4l In contrast, for symmetric 
noise the dynamics is unaffected, no matter whether cou- 
pling and phase asymmetries are present or absent. The 
only difference is that the incoherent state becomes un- 
stable for larger values of the critical parameters as one 
increases noise intensity. 



C. The effect of phase asymmetry 

For the case a 7^ 0, the inter-ensemble regions Dl and 
SI shrink as a increases, whereas the intra-ensemble syn- 
chronization region S2 expands, as shown in Fig. [71 This 
means that finite phase asymmetry reduces the probabil- 
ity of inter-ensemble synchronization (note reduced SI 
and Dl regions in Fig. [7]) and mostly allows only intra- 
ensemble synchronization of one or both of the ensembles. 
For a given set of parameters, on increasing a from 0, the 
following condition is satisfied 




m^D2 ^flf^i S1 



0.5 



B 



1.5 



FIG. 7: Theoretical B-a bifurcation diagram for (a) A'^' = 
1.2, A'-^^ = 1, (b) A'^) = 1.8, = 1.4 and Ac^ = 1,7 = 0.5. 
The line of *s represents the numerically obtained bifurcation 
boundary between the synchronized and incoherent states. 
Greatly reduced SI and Dl regions occur due to the presence 
of phase asymmetry. The discrepancy between numerical and 
analytic boundaries is discussed in the text. 



up to a critical value of a = aj given by 

'AAu; ± (8^2 - Aw2{i2 + b^))1 



sm 



4^ 



^ cos(2q;) + AAuj sin a - Auj^ > 



(13) 



where again, for a given set of parameters, there can only 
be one value of a that satisfies < a < 7r/2. Upon cross- 
ing aj, the condition is violated and the following 
condition is satisfied ^ cos(2a) -I- AAlu sin a — Aw^ < 0. 

As a result, when a > aj the inter-ensemble synchro- 
nization breaks down and the system enters into a state of 
intra-ensemble synchronization. Thus as one travels from 
SI (Dl) to S2 (D2) across aj the combined synchroniza- 
tion with single (double) frequency breaks between the 
ensembles and independent synchronization with single 
(double) frequency regime appears. Region S2, unlike 
region S in Fig. [I](a), embraces two states (i) synchro- 
nization in ensemble 1 with ensemble 2 incoherent and 
(ii) synchronization in ensemble 2 with ensemble 1 inco- 
herent, but does not distinguish between them. 

Further, there is a critical value of a = ar above which 
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FIG. 8: The coherence parameters r^^-* (grey) and r'^' (black), 
and the phase difference 5xp, plotted as functions of time, 
obtained from numerical simulations. Parameter values were: 
top panel, B = 0.7, a = 0.2; middle panel, B = 0.7, a — 7r/4 
(near Qc); and bottom panel, B — 0.7, a — 1.2 {a > etc), 
corresponding to regions D2 near the NS/S2 boundary and 
NS respectively (see Fig. Oa)). 



the collective oscillations disappear and the incoherent 
state becomes stabilized (see Figs. [7]and|S]). Thus by 
reducing the chances of occurrence of inter-ensemble syn- 
chronization and favoring intra-ensemble synchroniza- 
tion, phase asymmetry plays a crucial role in determin- 
ing the route to synchronization. Thus, for instance, in a 
particular problem, if one wants to have (avoid) the phe- 
nomenon of inter-ensemble clustering (region Dl) then it 
is obvious that phase asymmetry should be absent (finite, 
large). 

The discrepancy between the numerically and analyt- 
ically obtained boundaries in Fig. [7] is attributable to 
the influence of phase asymmetry. This affects region 
S2 which is large here (cf. Fig. Hfb) where both S2 and 
the discrepancy are smaller) and it changes the thresh- 
olds for A^") and Note that neither numerics nor 
analytics provides an exact result. The analytic bound- 
ary is obtained from the condition A^'^') > and refers 
to the limit of infinitely many oscillators. For numerics, 
the (asterisked) boundary is obtained from the condition 
j.(i.2) ^ Q j g^j-^j refers to a finite number of oscillators. 



D. Stability of the fully synchronized states in the 
limit 7 = 

In this subsection we focus on the noise-free case with a 
frequency distribution that has an infinitely sharp peak. 
In this case, the dynamics is reduced to that of two en- 
sembles of identical oscillators (the oscillators within the 
ensembles are identical while the ensembles themselves 
are non-identical). Now the Sip corresponding to intra- 



ensemble synchronization can be obtained from Eq. ([T]) 
as 

= sm~^[{Auj — j4sinQ:)/2i?cosa]. 

A linear stability analysis of Eq. ([T]) then gives N — 1 
degenerate eigenvalues for each ensemble, namely 

A± = - A^^'^^ cos a - S cos(±(5V^ -I- a) < 0, 

that characterize the stability of the intra-ensemble syn- 
chronized states of ensembles 1 and 2 respectively. In ad- 
dition, two eigenvalues Aq = and Ac — ~2B cos a cos Sip 
characterize the stability of inter-ensemble synchroniza- 
tion. Hence the transition between inter-ensemble and 
intra-ensemble synchronization states occurs at the fol- 
lowing bifurcation point 

Be — {Aoj — A sin a) /2 cos a. 

Note that 6ip varies from — 7r/2 to 7r/2 as a increases. 
For a < ah (not shown in figures), the stability condi- 
tion is satisfied by both ensembles and so intra-ensemble 
synchronization occurs in both ensembles. When a > ah 
the stability condition is violated by either one of the en- 
sembles and at that point a Hopf bifurcation occurs. As 
a consequence, intra-ensemble synchronization occurs in 
one of the ensembles. For the case ^(1) = A<-^\ when SiP 
varies from 7r/2 to for increasing a and when a > ah 
the stability condition is violated by the first ensemble. 
On the other hand, when Sip varies from to —n/2 with 
increasing a, the stability condition is violated by the 
second ensemble above ah- 



V. ROUTES TO SYNCHRONIZATION 

Given that the system possesses distinct synchroniza- 
tion regimes, it is of interest to investigate the routes 
it follows to synchronization. As one would expect, the 
route depends on the coupling, noise and phase asymme- 
tries. In particular, in the presence of coupling asymme- 
try, we have identified the following routes grouped 
into the two different cases a = and a ^ 0, and as- 
suming that we increase the inter-ensemble coupling pa- 
rameter B keeping all the other parameters fixed. When 
a = we find that there are at least three typical routes: 

1 . The oscillators in the ensembles pass from the syn- 
chronization regime D2 through Dl to the region 
SI. Thus when the ensembles are synchronized sep- 
arately, increasing B results in inter-ensemble clus- 
tering which then leads to inter-ensemble synchro- 
nization or phase locking between the ensembles. 
This route is represented by line TII of Fig. [IJb). 

2. There is also a possibility that when the ensem- 
bles are synchronized separately and when B is in- 
creased, the intra-ensemble synchronization be de- 
stroyed in one of the ensembles, which on further 
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increase of B, leads to phase-locking between the 
ensembles. Thus when the system is in region D2 
increasing B causes the system to pass through the 
region S2 to region SI. Inter-ensemble clustering 
does not occur in this route. 

3. If the ensembles arc initially not synchronized (that 
is in region NS) , then increasing B can cause phase- 
locking of the ensembles directly. Thus the system 
can pass directly from region NS to to SI. This 
route is characteristic of the case a = and cannot 
occur in the presence of phase asymmetry. 

In the presence of phase asymmetry, i.e. a ^ 0, the en- 
sembles can follow any of the following routes to synchro- 
nization 

1. The ensembles pass from region D2 through Dl to 
SI. This route is similar to route 1 that occurs 
for the case a = 0. Note that when A^^) = 

or Alu = only one entrainment frequency exists 
below 7c- and therefore this route does not occur 
for either cases (due to the non-occurrence of region 
Dl). 

2. The ensembles pass from region D2 through S2 to 

51. This route does not incorporate the state of 
inter-ensemble clustering. 

3. When the ensembles are synchronized separately 
(in region D2), increasing B causes the disruption 
of synchronization in one of the ensembles leading 
to synchronization in the other ensemble (region 
S2). Thus the ensembles pass from region D2 to 

52. This route is characteristic of phase asymmetry 
and cannot occur for the case a = 0. 

4. If the ensembles are not synchronized, increasing 

B will result in synchronization in cither one of the 
ensembles. Thus the ensembles pass from regions 
NS to S2 (unlike NS to SI in the absence of phase 
asymmetry). 

Knowledge of these routes to synchronization is obvi- 
ously important for the control of synchronization in real 
systems. 

VI. DISCUSSION AND CONCLUSIONS 

One might intuitively suggest that the synchronization 
phenomena induced by coupling and noise asymmetries 
could also be obtained by choosing a sufficiently large 
difference between the mean frequencies of the two en- 
sembles. However, the synchronization phenomena cor- 
responding to the Dl and S2 regions can only be ex- 



plained by introducing either coupling or noise asymme- 
tries. As an illustration let us consider the eigenvalue for 
the noise-free case without coupling and phase asymme- 
tries for Aw^ > B'^ 



A± = -7 + — ± i- VAw2 - _B2 - iu). 

For this case, the intra-ensemble synchronization takes 
place simultaneously in the two ensembles since the 
curves 7c+ and 7c- coincide when Re{X±) becomes pos- 
itive. Although there occur two Hopf bifurcations, they 
happen to be one and the same and hence one will not be 
able to explain the synchronization region S2. A similar 
problem occurs also with the Dl synchronization regime 
for Aw^ < B'^. Therefore we must conclude that the in- 
troduction of coupling/noise asymmetries are crucial to 
account for certain synchronization phenomena and can 
never be replaced by the introduction of large difference 
between the mean frequencies of the two ensembles. 

It is therefore essential to take account of possible 
asymmetry while attempting to model natural systems. 
Certain phenomena, like those discussed here, are at- 
tributable to asymmetries in the interactions. 

In this paper, we have investigated the role played by 
coupling, noise and phase asymmetries in two coupled 
phase oscillator ensembles. We have identified a global 
clustering phenomenon that may be characteristic of ei- 
ther the coupling or the noise asymmetry when the other 
is absent. Phase asymmetry reduces the likelihood of 
global clustering and also introduces new routes that are 
characteristic of itself. Thus phase asymmetry controls 
the routes to inter-ensemble synchronization. The phe- 
nomenon of inter-ensemble clustering that is characteris- 
tic of coupling asymmetry is found to occur even for sym- 
metrically coupled systems if noise asymmetry is present. 
Thus noise; asyninictry is found to complement the effect 
of coupling asymmetry when the latter is absent. 

We therefore conclude that, in modeling real sys- 
tems where synchronization arises, explicit consideration 
should be given to the effect of possible asymmetries in 
coupling, noise, and phase. 
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